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In an equilibrium system, the Kolmogorov-Sinai entropy, ^ks, equals the sum of the positive Lya- 
punov exponents, the exponential rates of divergence of infinitesimal perturbations. Kinetic theory 
may be used to calculate the Kolmogorov-Sinai entropy for dilute gases of many hard disks or spheres 
in equilibrium at low number density n. The density expansion of /iks is Ni'A\[nn + B + 0{n)\, 
where v is the single-particle collision frequency. Previous calculations of A were succesful. Calcu- 
lations of B, however, were unsatisfactory. In this paper, I show how the probability distribution 
of the stretching factor can be determined from a nonlinear differential equation by an iterative 
method. From this the Kolmogorov-Sinai entropy follows as the average of the logarithm of the 
stretching factor per unit time. I calculate approximate values of B and compare these to results 
from existing simulations. The agreement is good. 

PACS numbers: 05.45.Jn 



I. INTRODUCTION 

Chaotic properties of systems with many degrees of 
freedom, such as moving hard spheres or disks, have been 
studied frequently. Extensive simulation work has been 
carried out on their Lyapunov spectrum 0, |3j 13 1 ^^'^ ^'^^ 
low densities analytic calculations have been done for the 
largest Lyapunov exponent 0, IE IE > the Kolmogorov- 
Sinai entropy |E S, and for the smallest positive Lya- 
punov exponents ]E ^3 ■ 

The Kolmogorov-Sinai entropy describes the maximal 
rate at which the system produces information about its 
phase-space trajectory. In closed systems, it equals the 
sum of all positive Lyapunov exponents. 

In this paper, I consider a system consisting of N hard, 
spherical particles at small number density n, in d di- 
mensions [d = 2,3). I calculate the Kolmogorov-Sinai 
entropy in the low density approximation, where it is ex- 
pected to behave as 



/iKs = NDA [\n{na'^) + B + 0{na'^)] 



(1) 



where D is the average single-particle collision frequency, 
and a is the particle diameter. The constant A has 
been calculated by van Beijeren et al. but the re- 
sults found for B were unsatisfactory. I will show that 
the low-density approximation made in that calculation, 
which produces good results in the case of the Lorentz 
gas, which consists of uniformly convex scatterers, is too 
drastic in the case of a many-particle system. This has 
already been anticipated in reference 0| and preliminary 
estimates of the corrections were made by Dorfman in 
but the present calculations yield much more accu- 
rate values, which in principle could be improved even 
further. 



The paper is organized as follows: After a short intro- 
duction to Lyapunov exponents, the Kolmogorov-Sinai 
entropy will be introduced and its relation to the stretch- 
ing factor will be discussed in Sec. ^ followed by an 
explanation of the relevant dynamics of hard disks in 
Sec. mil In Sec. IIVI the stretching factor is calculated, 
and in Sec. IVlapproximations are introduced to the prob- 
ability distribution of the stretching factor. The paper 
finishes with a discussion of the results in Sec. IVII 



II. LYAPUNOV EXPONENTS 

Consider a system with an A/'-dimcnsional phase space 
F. At time t = 0, the system is at an initial point 7q 
in this space. It evolves with time according to 7(70, i). 
If the initial conditions are perturbed infinitesimally, by 
(57o, the system evolves along an infinitesimally different 
path 7 -|- ^7, which can be specified by 



(57(70, t) - M^,{t)-djo , 
in which the matrix {t) is defined by 



'70 



it) 



djQ 



(2) 



(3) 



The Lyapunov exponents are the possible average rates 
of growth of such perturbations, i.e., 



= lim In \fi.,{t)\ , 



(4) 
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where /ii(t) is the z-th eigenvalue of M^^, (t). If the system 
is ergodic, it comes arbitrarily close to any point in phase 
space for all initial conditions except for a set of measure 
zero. Therefore, the Lyapunov exponents are the same 
for almost all initial conditions. One may order the ex- 
ponents according to size, with Ai being the largest and 
AjV the smallest, as is the convention. For each exponent 
there is a corresponding eigenvector of Mjg{t). 
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For a classical system of A'' d-dimensional freely mov- 
ing hard spheres without internal degrees of freedom, the 
phase space and tangent space may be represented by the 
positions and velocities of all particles and their infinites- 
imal deviations, 



7^ = (rj,v,) , 



(5) 
(6) 



where i runs over all particles and Sji is the contribution 
of particle i to Sj. 



A. Kolmogorov-Sinai entropy and stretching factor 

In standard terminology, the stretching factor A(t) is 
defined as the factor by which the expanding part of tan- 
gent space stretches over a time t. This quantity can be 
used to calculate the Ruelle pressure as well as the sum 
of the positive Lyapunov exponents, which is equal to 
the Kol mog orov-Sinai entropy in systems without escape 
Ha. Il4l |. For long times, the stretching is dominated 
by the positive Lyapunov exponents, and one has for the 
Kolmogorov-Sinai entropy 




FIG. 1: Geometry of a collision of two particles in relative 
coordinates. The collision normal a is the unit vector pointing 
from the center of one particle to the center of the other. 



to collision i. In this paper, this includes the effects of 
the free flights of the colliding particles after the collision, 
and not those before. To obtain the Kolmogorov-Sinai 
entropy from this, one must calculate the distribution of 
single-collision stretching factors. 



'KS 



lim - In A(t) 

t^oo t 



(7) 



For long times, the stretching factor can be calculated 
from the total growth of an arbitrary volume element in 
dN dimensions. After a few multiples of the inverse of 
the smallest positive Lyapunov exponent, the dynamics 
project the volume element onto the expanding manifold 
and its subsequent growth is completely described by the 
stretching factor. 

For hard-sphere systems, where the collision times are 
defined exactly, the stretching factor can be written as 
the product of the stretching factors resulting from each 
of the different single collisions combined with the sub- 
sequent (or previous) free flights of the two particles in- 
volved. In this description, the effects of the free flights of 
the other particles are accounted for at the collisions in- 
volving those particles. On the right-hand side of Eq. iQ, 
the logarithm may be replaced by the sum of logarithms 
of these single-collision stretching factors. The resulting 
expression may be interpreted as a time average, which in 
ergodic systems may be replaced by an ensemble average. 
Hence, 



hKS = — (In A,; 



(8) 



At low densities, the single-particle collision frequency, 
9, is given by 



(9) 



The factor Nv/2 in Eq. © equals to the overall collision 
frequency. A; is the single-collision stretching factor due 



III. DYNAMICS OF HARD SPHERES IN 
PHASE SPACE AND TANGENT SPACE 

In order to calculate the single-particle stretching fac- 
tors, one must derive the dynamics of the system in tan- 
gent space from the dynamics in phase space. 

The evolution in phase space consists of an alternating 
sequence of free flights and collisions. During free flights 
the particles do not interact and the positions grow lin- 
early with the velocities. 



Vi{t) =r,(to) + (i-io)v.(to) 



(10) 
(11) 



At a collision, momentum is exchanged between the col- 
liding particles along the collision normal, a — {ri—rj)/a, 
as shown in Fig. ^ The other particles do not interact. 
Using primes to denote the coordinates in phase space 
after the collision, one finds 



I". = r^ 



0-(<T- Vy) , 



(12) 
(13) 



. jPll- 



where vij = Wi — Wj. From Eq. Q and Eos. jpn )l - l(T^ 
the dynamics in tangent space can be derived 

During free flight there is no interaction between the 
particles and the components of the tangent-space vector 
transform according to 



(5r' 



1 {t-h)i 

1 



Sr. I 
6wi 



(14) 



in which 1 is the d x d identity matrix. 
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A collision between particles i and j only changes the 
tangent-space vectors of the colliding particles As 
shown in Fig. ^ an infinitesimal difference in the po- 
sitions of the particles leads to an infinitesimal change 
in the collision normal and in the collision time. The 
v + (5v are exchanged along a + Sa according to Eq. (fT!^ . 
This leads to infinitesimal changes in both positions and 
velocities right after the collision. For convenience one 
may switch to relative and center of mass coordinates, 



Sri 



Sri - Srj,SRij = {Sri + 5rj)/2, = Svi - Svj, 



and SVij — ((5vj + S^j)/2. These transform as 



Sr'iA = Srin — 2S • Sr 



sn 



(5v^ = (5vy - 2S • - 2Q • Sr, 
(5V', = 5Y 



ij , 



in which S and Q are the d x d matrices 



(TCr , 



(15) 
(16) 
(17) 
(18) 



(19) 



Q ^ [(o- ■ Vjj) 1 + <7v,j-] ■ [{a ■ v,j) 1 - v.j-o-] 
a((7-Vy) 

Here the notation ab denotes the standard tensor prod- 
uct of vectors a and b. Note that Q transforms Srij vec- 
tors which are orthogonal to Vy into vectors which are 
orthogonal to v^^ . The vector Vy is a right zero eigen- 
vector of Q, and v^^ a left zero eigenvector. Note that 
these are d-dimensional vectors, not 2(i-dimensional. 



IV. STRETCHING FACTOR 

In order calculate the Kolmogorov-Sinai entropy from 
Eq. one must find the probability distribution of the 
single-particle stretching factor. In this section, I first 
derive the single-collision stretching factor as a function 
of the collision parameters and other parameters which 
contain information about the history of the system. In- 
formation about the history will be replaced by a pre- 
coUisional distribution function, which is the distribution 
function averaged over an ensemble of initial conditions. 



A. Projection 

The growth of a dA'^-dimensional volume element in ST 
can be monitored through its projection onto a subspace 
of ST with at least the same number of dimensions, as 
long as this projection space is not orthogonal to one of 
the dN leading eigenvectors of M. In the limit < ^ oo, 
the logarithm of the determinant of the transformation of 
the projection yields the same Kolmogorov-Sinai entropy 
as the logarithm of the stretching factor of the actual 
volume element. 

If {Sr^™^ , (5Vj-™'' ) are the eigenvectors belonging to the 
positive exponents, the eigenvectors which belong to 



their counterparts under conjugate pairing are equal to 
(Jrl™-*, — i5Vj^™-*). This means that eigenvectors which 
have no contributions along either Sri or Svi correspond 
to themselves under conjugate pairing. Such eigenvectors 
must therefore have Lyapunov exponents which are zero. 
The dA^-dimensional vectors whose components belong- 
ing to particle i are Sri respectively Svi are denoted by 
The spaces spanned by either Sr or Sv are not orthogonal 
to any eigenvectors which belong to nonzero Lyapunov 
exponents. Sr and Sw. In the system described here a 
convenient choice for the projection space may therefore 
be either of these spaces. Here, Sv is used for this pur- 
pose, because it does not change during free flights. 



B. Stretching factor of a single collision 

During a free flight, Sr grows with Sv. Denoting the 
perturbations in the position just after a collision with 
a superscript + and those just before the next collision 
with a superscript — , one may write 



Sr- — Srj 



(21) 



where r.^ is the free-flight time of particle i. Note that 
typically is of the order of 1/9. In previous calculations, 
it was usually assumed that right after a collision Sr^ 
and Svi were of the same order of magnitude Q. Under 
this assumption, the contribution from Srf to Sr~ may 
be neglected compared to TiS'Vi. Of course TiSvi will 
be comparable to Srf after the previous collision if Ti is 
short. However, this occurs only with a probability pro- 
portional to the density and therefore may be neglected 
in the average. 

The assumption that ^r^ and Swi just after a collision 
are of the same order of magnitude, however, is only true 
for d ~ I components of Srij , namely the ones normal to 
Vy. The remaining component of Srij, which is along 
Vij, and all components of SHij are, after a collision, 
larger by an order of r than the corresponding compo- 
nent of Svi, because Q, defined in Eq. H20|l . does not act 
on centre-of-mass perturbations, nor on perturbations of 
relative coordinates parallel to the relative velocity. In 
these directions, the components of Sv are of the same 
order of magnitude as before the collision, but the cor- 
responding components of Sr have grown linearly during 
the preceding free flights. I will show that this affects the 
Kolmogorov-Sinai entropy, even at low density. 

The determinant of the transformation of the dN- 
dimensional volume element projected onto Sw depends 
on Sr before the collision. Sr may be assumed to depend 
on (5v as 



Sr = tW ■ (5v 



(22) 



with f — Xjv the average free-flight time. W is pro- 
portional to the inverse of the radius of curvature ten- 
sor, which is often used to calculate Lyapunov exponents 
(see, for example, 23)- The matrix W can be split 
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up into d X d matrices between specific particles, \Nij. 
As particles collide and have free flights, W^j changes. 
The volume element projected onto before the colli- 
sion is mapped to a projection of a volume element after 
the collision. The determinant of this map depends on 
Wj, \Nj,\Nij, and Wji, where the second index is omitted 
if it is the same as the first. 

After the collision, W is changed. From now, a prime 
will be used to denote a quantity after a collision and 
the subsequent free flight(s). Quantities just after a col- 
lision, before the subsequent free flight are indicated with 
a superscript *. 

The matrix W after the collision, W*, can be found by 
using the dynamics and Eq. to express 6r* just after 
the collision in terms of Sv* , the collision matrices and 
yV. Let S and Q be the dN x dA^-dimensional matrices 
which perform the transformations of S and Q on the 
components of Sr and 6v along the colliding particles, as 
described in Eqs. H15|) - (|18|l and act as the unit operator 
on the components belonging to all other particles. The 
transformation can be written as 

6V = {J + S) ■ Sv + Q- Sr (23) 
= {I + S + fQ-W)-dv , (24) 

6r* = {J + S)■6r = f{J + S)■W■S^r . (25) 

Here, I is the dN x dN identity matrix. Note that 6v' 
is equal to 5v* , and that all components of Sr* and Sr' 
that do not belong to either of the colliding particles 
remain the same. This leads to an expression for Sr* as 
a function of S'v* . 



Sr* = fW* ■ Sv* , 



(26) 



where W* can be expressed in terms of W and the colli- 
sion matrices, as 



W* ^{I + S)-W-{I + S + fQ- 



(27) 



Using (2 + S) ^ = 2 + S, one may write this more con- 
veniently as 

W* = (I + S) ■ [W-^fil + S) ■ Q]-^ ■ (I + S) . 

(28) 

At low densities, two particles which collide can 
be assumed to be uncorrelated before the collision 
(StoBzahlansatz). This means that W^j ~ \Nji = 0, if 
i ^ j. After the collision there generally are nonzero 
elements in W*, . 

Wi also changes during the subsequent free flight. Let 
Ti be the free-flight time of particle i after the collision. 
Then, after the free flight. 



1 Wf 



ll/Tfe 



if k = i,j 
if k^i,j 



(29) 



Note that the change in Srk with k ^ i,j due to free 
flights was already taken into account at the collisions 



of particle k. The matrix W as it is calculated here de- 
scribes the connections between Sri just before the next 
collision of particle i and Svj just before the next collision 
of particle j. 

The matrix (2+S) ■ Q is nonnegative deflnite and sym- 
metric. If W is positive definite, so is its inverse. This 
means that W"^ -I- f(X -I- 5) • Q is positive deflnite, as 
is its inverse. The coordinate reflection X + 5 is unitary, 
and therefore the eigenvalues of {1 + S) ■ W* ■ {1 + S) 
are the same as those of W*. Therefore, W* is positive 
definite. By similar reasoning, a symmetric W is mapped 
onto a symmetric W* . Equation (|29|l also maps nonnega- 
tive deflnite matrices onto positive definite matrices, and 
symmetric ones onto symmetric ones. As, without loss 
of generality, any initial conditions for W may be cho- 
sen, it is possible to choose them such that W is positive 
definite and symmetric. This can be done, for example, 
by choosing the initial W to be diagonal, with elements 
equal to I. 

The stretching factor due to one collision can be cal- 
culated from the determinant of the transformation of 
the projection onto the perturbations of the relative and 
centre-of-mass velocities, Svij and SYij. From Eq. H24|l . 
one finds that this is the determinant oiT+f{T+S)-Q-W. 

In the low density limit, and with \Nij — 0, this is 
found to be equal to 



A = w±± 



2VT 



COS 



d-3 , 



(30) 



Here, 9 is the angle between a and v, cos 9 ~ ct • v, 
and w±±_ is equal to the determinant of the part of 
(Wi -I- Wj)/2 between vectors that are orthogonal to v.^ 
before the collision. For d — 2, 



(31) 



This expression replaces the factor Pr-)_, where r+ = 
(Ti + Tj)/2, in previous calculations of the single-particle 
stretching factor In d dimensions, must be re- 
placed by a set of d — 1 vectors orthogonal to v^j . The 
inner products are then replaced by the determinant of 
the (d — 1) X (c? — 1) matrix with elements given by the 
innerproducts of (W^ +\Nj)/2 between those vectors. 



V. DISTRIBUTION FUNCTIONS 

In order to calculate the Kolmogorov-Sinai entropy 
from Eq. ((SJ, one needs the distribution function of 
the single-collision stretching factor, as described by 
Eqs. H3U|) and H31|l . This may be derived from the joint 
distribution function of the collision parameters r^, Tj, v^, 
Vj, 9, and the elements of and W^. In the low-density 
approximation, the collision parameters are distributed 
according to the equilibrium solutions of the Boltzmann 
equation. 
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The distribution of the particle velocities is the 
Maxwell distribution, 

'/'m(v.) ^ (^^) ' exp (^-i/3m|v|2) . (32) 

The collision rate for collisions at angle Q and with outgo- 
ing velocities and Vj is proportional to the differential 
cross section times the relative velocity, sin''"^ OcosOvij. 
The normalized probability distribution of the collision 
parameters is thus equal to 

p{vi , Vj , 9) dwd^j dO ^ sin''-^ 9 cos 9 1^' " 

X (l)Miyi)<l>Miyj)d-Vid-Vjd9 . 

(33) 

The free-flight times of the particles are distributed ex- 
ponentially, with the collision frequency i/(v) depending 
on the velocity of the particle, according to 

Pr{T\Vi)dT = v{vi) exp[-i^(vi)T]dT , (34) 

where v{vi) is the velocity dependent collision frequency. 
This can can be found by integrating the collision rate 
over the outgoing velocity of the other particle and over 
the collision normal, i. e. 

ly(vi) — j dw jd9 na'^^'^ •Aa!^^'^ cos 9\vi ~ \j\4>M{"Vj) ■ 

(35) 

The distribution of can be found from the re- 

quirement that the distribution of elements of as a 
function of v.; is not changed by collisions. This yields 
a complicated nonlinear differential equation for the dis- 
tribution of the elements of and Wj . It involves the 
distribution of angles between the relative velocities of 
subsequent collisions of a particle with velocity as well 
as the velocity dependence of the collision frequency. The 
latter is only known numerically. With the implied inclu- 
sion of the collision parameters and coordinates of both 
particles in the distribution functions implied in the in- 
tegral, the equation can be written as 

p\W) = j dWp{W)6{W'{W) - W) . (36) 

The solution to this equation can be approximated using 
an iterative approach. 

A. Approximation of the distribution of w±± 

Rather than solving Eq. 1)36(1 exactly, which is not 
feasable, Eq. may be used to iterate the distribu- 
tion. One may start with an initial distribution, p{W), 
and calculate the distribution after one collision, p'(yV). 
After every iteration, the distribution more closely resem- 
bles the true solution of the integral equation. However, 



even with a simple initial distribution, such iterations will 
quickly produce distribution functions which can only be 
calculated numerically. In principle, the equation could 
be such that there is no convergence at all, but the nature 
of the physical problem dictates that the distribution of 
W converges after many collisions. 

In this section, an alternative iterative approach is used 
to find an approximate distribution function. We start 
with a simple distribution with one parameter, which ap- 
proximates the average trace element. The parameter is 
chosen in such a way that the average of the trace of W 
remains the same after a collision and free flights. The 
nonlinear terms in the equation for the distribution of the 
elements of W, which have zero average, can be ignored 
at first, due to the choice of initial distribution. The in- 
tegral equation is iterated a second time to include some 
of these terms. Subsequently, the size of the remaining 
corrections after more iterations is estimated. 



1. The trace ofW 

In principle, it would be possible to use the determi- 
nant or some other scalar function of W, instead of the 
trace. It is however much easier to write down the map 
of the trace of W onto the trace of W* in Eq. (j2Hl than 
it is to write down a map of the determinant during free 
fiights. Also, only the distribution of diagonal elements 
of yy is actually needed. Under unitary coordinate trans- 
formations, such as rotations and reflections, the trace of 
a matrix is conserved. Using a parameter is not really 
necessary, however, it greatly improves the convergence 
toward the solution of the equation. 

From Eq. I(29|l . the trace of W' can be found to satisfy 

Tr(W) = Tr(W*) + dvT, + dvT, . (37) 

Let the dTV-dimensional basis vectors in which the ma- 
trices are expressed be numbered 1 through dN . Let the 
first dA^-dimensional basis vector, ei, be defined as v.^ 
in the relative coordinates, and the second (in d dimen- 
sions the second through d-th) , £2 , as v^j _l in the relative 
coordinates. The remaining basis vectors may be chosen 
in any arbitraty way, as long as they are orthogonal to 
eachother and of unit length. 

Let Ai represent the dN x dN matrix A with all rows 
and columns removed except for those with indices speci- 
fied by the list where I may be any list of indices. Sim- 
ilarly, let A(^i-) be the matrix A with all rows and columns 
removed with indices belonging to the list I. Specifically, 
A(k) and A(k2) represent the matrix A from which the 
rows and columns belonging to, respectively, index k and 
both index k and index 2 are removed. 

From Eq. H28|l . the trace of W* can be found as a 
function of the elements of W. As the trace is conserved 
under the coordinate reflection {X + S) -W ■ [I + S), one 
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finds 



Tr(>V*) = 



dct(W-^ +f(J + 5) ■ 
det(>V-i +f(X + 5) • Q) 



(38) 



Here, the sum is over an ortlionormal basis of dN unit 
vectors. In the low-density Umit the mean free time 
becomes large, and only terms in which the numerator 
contains the same power of r as the denominator can 
contribute. The product of the nonzero eigenvalues of 
{1 + S) ■ Q can be divided out, leaving only the deter- 
minant of the remaining part of W, between vectors on 
which Q does not work. As only {T + S) ■ Q in Eq. 138|) 
contains the collision normal, the trace of W* docs, in 
the limit of vanashing density, not depend on 6, but only 



on Vij and the elements of W. 



The trace of W* can be rewritten as the sum over 
fractions of subdeterminants, 



Tr(>V*) = J2 

k^2 



det(ng)) 
det(W(-i) 



(39) 



Further, after writing the inverse of W^^ and by working 
out the determinant of W^, by induction over the number 
of indices occurring in i 



detW(-)i 



detW^ detW"! . 
From Eqs. (gUl and one finds, for d = 2, 
det(Wfe2) 



Tr(>V*) = 5] 



detWa 



k^2 

Tr(W) - 



1 



(40) 



(41) 



(42) 



In the d-dimensional case, £2 is replaced by d — 1 vectors, 
and Eq. (|42() becomes somewhat more complicated. 



The change in off-diagonal elements at a collision can 
be found from a derivation similar to that for the trace 
in Eq. (O, 



el ■ W*t: = ep • W • e 



(ep->V-e2)(e2->V-£g) 
£2 • W • e2 



(43) 



The expression for d = 3 is similar. If p or (7 is equal to 2, 
the off-diagonal element vanishes. Off-diagonal elements 
between different particles are not affected by free flights, 
as is apparent from Eq. if^ . 

The collisions are most conveniently expressed in the 
basis which consists of v^j and the d — 1 vectors Vij_L 
orthogonal to it. I therefore also express each \Nij in this 
basis. 



2. Iterative approach 

Assume that just before a collision the are equal 
to their averages and W^j all zero. If the distribution of 
angles between the relative velocities of two consecutive 
collisions is (nearly) isotropic, the two average diagonal 
elements are (approximately) equal, so that 



(44) 



where 5ij is the Kronccker delta. The initial distribution 
used in the iteration process is a product of Dirac delta 
functions at the average value w for the diagonal elements 
and at zero for the off-diagonal elements. In a similar 
way, an exponential distribution function can be used, 
with average to test the sensitivity to the width of the 
distribution. 

Using Eqs. (|28|l and (|29|l . one finds that after the col- 
lision and free flight, in the basis consisting of and 
the d—1 vectors orthogonal to it, the values of Wfej have 
changed according to 



kl 



[W + VTk) 



{\w + VTu)ld~ 




w 1 Ski 



if k^l^iyk = l^j 



if (fc,0 = (*,j)v(fc,0 = (j,*) , 



if k^iJVl^iJ 



(45) 



where Id-i denotes the {d— l)-dimensional identity ma- 
trix. This equation implies a distribution for the elements 
of expressed in the basis belonging to the next colli- 
sion, which consists of v^^ and the d—1 vectors orthogonal 
to it, v.^y^. The new distribution of the matrix elements 



is the distribution of \Nki in the coordinates of the next 
collision, TZi ■ WJ,; • TZf , where TZi is the rotation matrix 
associated with the rotation from the coordinate system 
using the post-coUisional relative velocity of a collision 
to the system using the pre-collisional relative velocity of 



7 



the next collision of the same particle. In two dimensions, 
this matrix is characterized by the angle between the 
relative velocities at the two collisions, 



7^. 



COS0 

— sin ( 



sm(; 
cos ( 



(46) 



The distribution of this angle depends on the velocity of 
the particle between the two collisions. In three dimen- 
sions, the angle 0i is to be replaced by two angles. 

From this approximation a distribution function of 
w±± can be found, which depends on w. At the next 
collision, it is, for d = 2, equal to the distribution of 



v^.(w^ + w^).v^ 



(47) 



= -(0, 1) • (7^. • . 7^f + ■ w' • nj) ■ (o, i) 



1 (cos^ I 

4 



cos 



(48) 
(49) 



The distribution of w±± can be approximated by the 
distribution of the right-hand side of the equation. For 
d = 3, the rotation matrix is more complicated, and so 
Eq. (|49|l becomes more complicated. The resulting ex- 
pressions are not reproduced here. 

From Eq. H45|l one can find the difference between the 
average traces of W" and W. One finds, 



(Tr(>V')) - (Tr(>V)) = 2d - {d - l)i 



(50) 



where the notation (.) denotes the ensemble average. As 
the average trace must not change, one finds for w the 



approximation 



r,(0) 



2d 
d- 1 



if d = 2 , 
if d = 3 . 



(51) 



Note that this result for id would be the same if an ini- 
tial exponential distribution were used for the diagonal 
elements. The resulting distribution function for w±±, 
however is different in that case. The distribution func- 
tion implied by Eq. (|49|l for w±± at the next collision 
can be used to estimate the Kolmogorov-Sinai entropy 
from Eqs. ((SJl and (|3U|I . The involvement of these expres- 
sions was already anticipated in pj] | by Dorfman, who 
predicted extra contributions to B of In 4 and In 3 based 
on summations. 

The approximation so far is fairly crude. The nonlin- 
earity of Eq. (|28|l has been partially neglected by using 
the averages of the off-diagonal elements. In the second 
term in the calculation of the trace in Eq. (|42|l only the 
block diagonal terms, those of the form v^jx • W| • v^j., 
are involved. In reality, since W is symmetric, the terms 
involving off-diagonal elements will also produce negative 
contributions to the average of the trace in Eq. H42I) . 

A better approximation of the average value can be 
found by iterating the equation for the distribution a sec- 
ond time. The distribution of TZi ■ \N'f.{Rj can be used 
to calculate the trace of W". The colliding particles are 
uncorrelated before the collision, but not independent of 
the particles they encountered before. These particles, 
which are not directly involved in the collision, now con- 
tribute to the change in the trace, through the second 
term on the right-hand side of Eq. (|42|) . I find that for 
d = 2 the trace in the low density approximation satisfies 



(Tr(>V")> - (Tr(W)) = -3u; + 

/ \2v'^{Tf + ATiTj + t|) + l{)v{Ti + Tj)w + 8w'^ — 2'Wv{Ti COS^ 4>j + Tj COS^ — W^(C0S^ + COS^ 



■ 4w 



(52) 



This yields a result for w that is significantly different 
from Eq. if^ . 

More iterations would produce more terms and will 
further reduce the value found for iD, converging to the 
exact result. A similar but far more complicated expres- 
sion can be found for d = 3 from Eqs. (|42|) . (|45|) . and the 
general form of the three-dimensional rotation matrix. 
The results would improve if the distribution were iter- 
ated repeatedly, but this would produce expressions of 
complexity increasing exponentially with the number of 
iterations. One more iteration would add four angles and 
four free-flight times to the expression in Eq. (|52|l . Af- 
ter the second iteration, the expression for ui" j_ is quite 
complicated and contains twelve correlated variables, six 



rotation angles and six free flights. Also, the second it- 
eration already produces a reasonable result. I there- 
fore continue using the distribution of ui^j^, but with the 
value of w found from Eq. (|52|) . 

The integrations over the distributions of (/i^, 0j, r^, and 
Tj can be done numerically. The change in the trace is 
zero for for w equal to 

(i)_r3.009 ifd = 2, 

"\2.107 if d = 3. ^^^^ 

The subscript index is introduced to indicate the weight 
given to the off-diagonal terms. 
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s q r 




FIG. 2: A diagrammatic representation of a series of col- 
lisions. At the collision between particles q and r, the off- 
diagonal elements between q and particle s contribute to the 
change in the trace. These elements date from the collision 
between particles x and y, in the common history of q and 
s. At some point in the past, there were elements between a 
particle i in the history of s and a particle k in the history of 
q when j, another particle in the history of s, collided with 
i. After the collision there were elements between j and fc, 
which, through more collisions, eventually lead to elements 
between s and q. The size of the elements between s and q 
can be estimated using Eqs. 1591 and Ifcim . 

3. Off-diagonal elements from earlier collisions 

If the contributions to the trace from the off-diagonal 
elements involving the other particles from the previous 
collisions, through the second term on the right-hand side 
of Eq. 142|1 . are ignored, the result is changed signifi- 
cantly. In this case, 

a)^r3.408 ifd = 2, 

° [2.639 if rf = 3 . ^ ' 



At a collision between i and j, the off-diagonal elements 
between particles i and k produce significant changes to 
the diagonal elements of WJ,,;,. It is therefore expected 
that contributions from particles involved in collisions 
before the previous collision will also be significant. Also, 
if the other particle from the previous collision of particle 
i has collided since, this has an effect on W^. 



At a collision between particles q and r, the change in 
the trace of W, calculated in Eq. (|42|l . is affected by the 
elements of W between q and r, and other particles. Af- 
ter a collision between two particles, nonzero off-diagonal 
elements exist between these particles. After a collision 
between i and j, off-diagonal elements between particles 
i and k generate off-diagonal elements between j and fc, 
due to the exchange between the Jv^ and 5vj . If nonzero 
off-diagonal elements exist between i and k as well as j 
and I before the collision, after the collision nonzero ele- 
ments will exist between k and I. A diagrammatic repre- 
sentation of the collision sequence is shown in Fig. |2] 



In order to estimate how much such terms contribute 
to the change in the trace at a collision involving parti- 
cle k, the typical magnitude of the off-diagonal elements 
at a collision must be investigated. One may estimate 
the typical changes in the off-diagonal blocks WjA: and 
Wj-fe at a collision between i and j, by estimating the 
changes in the trace of the off-diagonal blocks. The typi- 
cal size of the off-diagonal elements can be characterized 
by the trace of the off-diagonal block Wife. The diago- 
nal elements of the off-diagonal blocks can be found from 
Eq. Using the fact that Vslij before the coUision is 

zero, one finds, in two dimensions, for the elements of W 
between 5vi and i5vfe 



Vy • W-fe • e = • Wife • e - 



(vjj ■ W, ■ v.j-±)[v,j-± ■ (W,fe - Wj-fe) 



v.,-L • W^fe 



(W,fc - W,fe) 



,i . (W, - W,) • v,ju][vy-± • (W,fe - W,fe) 



2v,: 



(W,+Wj).v„i 



(55) 
(56) 



Here, e can be any vector in two dimensions. If Wife has 
nonzero elements, then W^-fe does not, since the particles 
i and j were uncorrelated before the collision. If W^fe has 
nonzero elements after the collision, both W^^, and W^^, 



have nonzero elements. 

From Eq. (|55|l the traces after the collision may be 
found. 
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Tr(W:, 



Tr(W' 



jk) 



• Wife • 



{vij ■ \Ni ■ ■Vijj_){vij± ■ \Nik ■ Vij) 



1 - 



^ aaT — I \A/ \ — : oi^u-L ■ VVife • Vy-j 



r 



(57) 
(58) 



It is fair to assume that the off-diagonal elements of W^fc 
tend to be smaller than the diagonal elements. Also, the 
diagonal elements of W — Wj are, typically, much smaller 
than the diagonal elements of + \Nj . The terms with a 
quotient of these can therefore be neglected in this rough 
estimation. Further, the d diagonal elements of W^fc are 
typically of the same size. These approximations leave 
us with the general expressions 



Tr(W:,) « ^Tr(W,,) , 



Tr(W;., 



2d 
d-1 
2d 



Tr(W,fc) 



(59) 
(60) 



In addition, if both W^fe and Wj; have nonzero ele- 
ments, WJ.; also has nonzero elements, which are due 
to the second term on the right-hand side of Eq. (|43|l . 
The off-diagonal elements generated in this way are small 
compared to the elements generated from Eqs. (|59ll and 
PHl . In fact, they are smaller than the terms neglected 
from Eqs. and because they contain products 
of the off-diagonal elements which are small compared 
to the diagonal elements. As they appear quadratically 
in the change in the trace [see Eq. (|42l) ]. they may be 
neglected, despite their quadratically larger number. 

From Eqs. (|59|l and H60|l an estimate can be made of 
the contributions to the change in the trace in Eq. H42|) 
from collisions before the previous collision, compared to 
the contributions from just the previous collisions. The 
ratio between the total contributions from off-diagonal 
elements to the change in the trace and the contributions 
from just the previous collisions is denoted by a. 

In every collision in the history, off-diagonal elements 
are created between the two colliding particles and ex- 
isting elements are reduced in magnitude and passed on 
according to Eqs. H59|l and H6(J|I . In order to estimate the 
consequences of W^s at a collision between q and r (see 
Fig. . 121) , one has to find the path through which informa- 
tion was passed on from the collision between particles 
X and y to the present collision between q and r as well 
as the path from the collision between x and y to the 
particle s at the time of the collision between q and r, 
following a sequence of collisions, through which the off- 
diagonal element between particles q and s is affected. 
To this path belongs an approximate reduction of the 
size of the off-diagonal element, a product of factors of 



{d+l)/ {2d) or (d—l)/ {2d) for each collision in the paths. 
If the path continues with the same particle, there is a 
factor of {d+ l)/{2d). If it switches to the other particle, 
the factor is {d- l)/(2d). 

Every different product with the same number of fac- 
tors follows a different path of that length, and hence 
belongs to a different present particle. The product of 
the two factors of two paths starting from x and y gives 
the order of magnitude of the off-diagonal element be- 
tween i and k particles. The square of this factor then 
gives the relative size of the contribution to the trace at 
the collision between i and j. If a collision between two 
particles is now p collisions ago, then, on average, the 
other part of the future of that collision has also had p 
collisions. Summing over all the different paths of length 
p, one finds that the relative contribution from collisions 
that occurred p collisions before the previous collision can 
be approximated as 



d^ + 1 
2d2 



2p 



Summing over all p gives the estimate 
p 



(3d2 + l)(d2-l) 



(61) 



(62) 



The contributions from the previous collisions in 
Eqs. (|^ and (|^ can be multiplied by a, to find an 
estimate for the total contribution of all particles with 
which i and j have a common history. This is an admit- 
tedly crude estimate, yet should give better results than 
just neglecting the history before the previous collision. 
The terms in Eq. (|52|l that are due to the off-diagonal 
block between the colliding particles and other particles 
may be multiplied by a. With this correction it is found 
that 



w, 



(1) 



2.929 
1.947 



ff d: 

ff d- 



2 , 

3 . 



(63) 



The distribution function of w±± contains an uncer- 
tainty in its width, which affects the results of the cal- 
culation. When starting from the average, with every 
next iteration of the equation for the distribution func- 
tion, the distribution becomes wider. w±± looks like a 
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sum of several weighted free-flight times for each particle. 
If one starts from exponentially distributed diagonal ele- 
ments, rather than simply the averages, the distribution 
becomes narrower with every iteration. By starting from 
an exponential distribution, one may estimate the conse- 
quences of the width of the distribution of the elements. 
With an initial exponential distribution, one finds 



In earlier calculations [g, A was calculated accurately. 
The results for B, however, are unsatisfactory. The val- 
ues of A are easily found from Eqs. ijHl, PHjl . and the 
dependence of the collision frequency on n, 



d-1 



(66) 



r,(l)e _ 



2.426 



if d = 2 , 



1.676 if d = 3 



(64) 



If w±± is taken equal to t+, the results for B of Q are 
reproduced. 



By substituting the distribution fmiction induced by 
one iteration, together with the average, into Eq. ((SJ one 
can now estimate the Kolmogorov-Sinai entropy. 



hKS = ~2~ \ 



d-1 



„d-3 , 



(67) 



VI. RESULTS AND DISCUSSION 

As the free-flight times are inversely proportional to 
the density, w^^f'^~^ will be inversely proportional to 
^jj^jg leads to a general form for the Kolmogorov- 
Sinai entropy. 



/iKs = NvA[- \n{na'^) + B] . 



(65) 



This yields 



B 



0.209 if d = 2 
-0.583 if d = 3 



(68) 



From molecular-dynamics simulations the 
Kolmogorov-Sinai entropy has been calculated 0, 0| . It 
is found that 



(0.499 ± 0.001)iVi^ (- Inna'' -I- 1.366 ± 0.005) 
(1.02 ± Qm)Nv (- Inna'* + 0.29 ± O.Ol) 

I 



if d = 2 
if d = 3 



(69) 



In the calculation presented here, Eq. (|67|l has to be 
amended, to become 



Nv 



^KS 



In 



d-l 



„d-3 , 



(70) 



From Eqs. (|49|l and (|53|) . one finds, after numerical inte- 
gration, that 



(1) 



1.592 if d = 2 , 
0.476 if d = i . 



(71) 



If the contributions from the off-diagonal elements in 
Eq. (|42(l are increased by the estimate of the remaining 
terms, through a factor of a, the results change to 



5(1) 



1.572 \i d^2 
0.427 if d = 3 



(72) 



This more closely reproduces the simulation results 
shown in Eq. (|^ . 

After every extra iteration in the calculation, the dis- 
tribution becomes wider and therefore the average of the 



logarithm of w_\_^_ becomes smaller compared to the loga- 
rithm of the average. Due to this, cutting off the process 
after two iterations produces a result for the Kolmogorov- 
Sinai entropy which is too high. Note that also a wider 
spread of the off-diagonal elements leads to larger con- 
tributions from the off-diagonal terms in Eq. H42|) . and 
therefore to a smaller iD, which yields a smaller value for 
B. Equation ifTTjl threfore gives an upper bound for B. 

By starting from a wider distribution of diagonal ele- 
ments instead of a product of Dirac delta functions, an 
estimate can be made of the effects of the width of the 
distribution. From the results of an exponential initial 
distribution, Eq. an estimation is found. 



B(l)e 



1.370 if d = 2 , 
0.273 if d = 3 . 



(73) 



From these two estimated bounds, a final estimate of 
B may be made, including error bounds, 



B = 



1.47±0.11 
0.35 ±0.08 



\i d = 2 
if d = 3 



(74) 



The errors could be reduced by using distribution func- 
tions for w^i^ that have been iterated a larger number 
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of times. The values of the Kolmogorov-Sinai entropy 
found in the molecular dynamics simulations, Eq. (|69|1 . 
are well within the error bounds of Eq. (j74|l . 

VII. CONCLUSIONS 

In this paper, an estimation using kinetic theory is pre- 
sented of the Kolmogorov-Sinai entropy of dilute hard 
sphere gases in equilibrium. Kinetic theory has been ap- 
plied before to calculate chaotic properties 0, 0, lrn |. 
such as the Lyapunov spectrum of the high dimen- 
sional Lorentz gas systems with escape, the 
Kolmogorov-Sinai entropy has also been connected to 
transport coefficients [Ta Il9ll2(ll2ll| . but in the system 
investigated here, it is just equal to the sum of the pos- 
itive Lyapunov exponents. It is known that the leading 
orders are of the form NDA{— \nn + B)^. 

An nonlinear integral equation was derived for the joint 
distribution function of the elements of the inverse of the 
radius of curvature tensor. This equation was approxi- 
mately solved by the use of an iterative method. B was 
estimated from the solution in a satisfactory way, with re- 
sults which are consistent with simulation results. It was 
found that, B = 1.47±0.11 for d = 2 and B = 0.35±0.08 
for d = 3. The approximations made are systematic, and 
the results can be further improved by performing more 



iterations of the distribution. The values for B found in 
the present calculation are in good agreement with the 
results from molecular-dynamics simulations Also, 
an upper bound was found for B, that is, B < 1.592 for 
d = 2, and B < 0.476 for d = 3. 

The smaller Lyapunov exponents of this system are 
proportional to P. The ones which are not due to 
Goldstone modes contribute significantly to B in the 
Kolmogorov-Sinai entropy. The calculation of B pre- 
sented here shows effects which affect the behavior of 
these Lyapunov exponents to the leading order [23| . 

It should be noted that effects such as the ones de- 
scribed here do not affect the Lyapunov spectrum of 
the high-dimensional Lorentz gas, calculated in reference 
|l2j| . because the scatterers in that system are uniformly 
convex. However, they are generic for the Lyapunov spec- 
tra of systems consisting of many particles. 
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